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Abstract 

We consider a three-state model comprising tumor cells, effector cells and tumor detecting cells 
under the influence of noises. It is demonstrated that inevitable stochastic forces existing in all three 
cell species are able to suppress tumor cell growth completely. Whereas the deterministic model 
does not reveal a stable tumor-free state, the auto-correlated noise combined with cross-correlation 
functions can either lead to tumor dormant states, tumor progression as well as to an elimination 
of tumor cells. The auto-correlation function exhibits a finite correlation time r while the cross- 
correlation functions shows a white noise behavior. The evolution of each of the three kinds of 
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O cells leads to a multiplicative noise coupling. The model is investigated by means of a multivariate 



Fokker-Planck equation for small r. The different behavior of the system is above all determined 
'""' by the variation of the correlation time and the strength of the cross-correlation between tumor 

^ and tumor detecting cells. The theoretical model is based on a biological background discussed in 

Tjj- detail and the results are tested using realistic parameters from experimental observations. 
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I. INTRODUCTION 

Tumor growth has become an important issue in medicine, biology and physics. The un- 
derstanding of cancer growth mechanisms is necessary to develop relevant strategies against 
the disease. In the past, deterministic models have been proposed for interacting tumor and 
immune cells which are investigated by performing stability and bifurcation analysis [lH3]. 
Moreover, a deterministic mathematical model with strong relation to experimental data is 
presented in [4j. As a new aspect the delay time between the detection of tumor cells by 
the immune system and the arrival of activated killer cells at the tumor site was taken into 
account in [S]. All these mathematical models can be considered as two state models of 
predator-prey-type. In general, such models can show interesting behavior as demonstrated 
in many examples in [6]. Recently Ref. [7] has discussed the effect of deterministically 
imposed transitions in reaction and population systems on the rates of rare events such as a 
crossing-over to population extinction. Another approach was chosen in [8] where the early 
stages of tumor growth was investigated. More precise, the geometrical aspect of contour 
instabilities was related to cell-cell interactions. Likewise the role of noisy influences can be 
regarded. As a result the stochastic forces may change the dynamics, in particular it was 
shown that the evolutionary dynamics is altered in case demographic noises are included in 
a deterministic model of interacting players |9]. As well, intrinsic stochasticity was consid- 
ered in [To] applied to the Lotka-Volterra model with special emphasis on the elimination 
of species. In addition, the extinction of stochastic populations caused by intrinsic noise 
was analyzed in fH]. Regarding tumor evolution one often refers to a logistic growth model 
which offers relevant results in spite of its simplicity p^. In the present paper we also use 
as the basic model the logistic equation for the deterministic cancer cell growth dynam- 
ics, see Eq. ([I]). Recently a generalized logistic equation was studied by supplementing the 
birth rate by a Markovian dichotomic noise [13]. Another essential point is that the tumor 
genesis is often accompanied by an abnormal proliferative activity of human tissue. In [TJ] 
the authors have reported on a mathematical model which covers the growth properties in 
terms of a variable renewal rate of cell populations in colon crypts. A further class of models 
is related to a single population where only the tumor cells are considered as the relevant 
variable. Here the deterministic equation are subjected to additional random forces which 
allows an analysis in terms of the related Fokker-Planck equation. Models for for white noise 



[15] and for colored noise [12] has been predicted. The latter one contains tumor-immune 
interactions in an implicit manner. Later a modified model was investigated by introducing 
a bounded noise which mimics the reduction of the tumor size due to a possible immune 
response [16]. Therefore, the random nature is also attributed to the immune system. 
This idea plays likewise a significant role in the present approach. Different to former works 
[12] [T6] the tumor-immune interplay is now incorporated explicitly. However the main point 
is that we demonstrate tumor-immune cell reactions can be induced by stochastic forces. 
To be more specific our model describes the time evolution of three different cell types: (i) 
tumor cells the density of which is denoted by X{t), (ii) effector cells with density Y{t) and 
(iii) tumor detecting cells with density Z(t). Whereas the last kind of cells is only able to 
recognize tumor cells but not to kill them. The effector cells have the ability to eliminate 
tumor cells. The deterministic model introduced in Sec. Ill describes the mutual interaction 
between the three species. However this model offers no stable tumor-free state. Due to the 
inclusion of inevitable randomness the growth and death rates of the immune and tumor 
cells, respectively, are altered immediately. Toward a more realistic description we allow also 
the occurrence of cross-correlations between the noise acting on the tumor cells and that one 
acting on the detecting cells. The resulting set of stochastic equations with multiplicative 
noise can be transfered to the related Fokker-Planck equation. By variation of the strength 
of the cross correlation and the finite the correlation time the system tends to different stable 
states which differ from those of the deterministic system. Especially we show that the noisy 
system exhibit the complete suppression of the tumor. The paper is organized as follows. 
First, we present in Sec. II some biological ideas concerning the tumor-immune interaction 
The mathematical model is developed and discussed in Sec. III. Due to the inclusion of ran- 
domness the related Fokker-Planck equation is introduced in Sec. IV. Afterwards we present 
our results in Sec. V before we finish with some conclusions in Sec. VI. 

II. BIOLOGICAL MOTIVATION 

Before we present in the forthcoming section the mathematical model let us summarize 
some biological mechanisms concerning the interaction between the tumor and the immune 
system. In particular, this section is focused on the main underlying ideas which are nec- 
essary for the understanding of our presented model. Introduced in the early 1900's [17] 



and again suggested in the middle of the 20th century [181 HH] there is the hypothesis that 
the immune system is able to detect and to eliminate nascent transformed cells. During 
the last decade the concept of the immune surveillance of the body played a significant role 
in tumor elimination, too. The investigations are supported by experimental results verify- 
ing the immune surveillance hypothesis |20](and Refs. therein). Furthermore, the immune 
surveillance concept was modified and is now known as 'immunoediting' which reflects the 
dual role of the immune response during the early stages of cancer growth [2T1 - I23] . The term 
immunoediting means both the ability of the immune system to destroy the tumor cells and 
a possible sculpting of the cancer cells. As the result all cells with a low immunogenicity will 
survive and begin to proliferate. This escape of the tumor from the control of the immune 
system can be regarded as a special feature of tumor growth [21] . As the consequence of the 
transformation of normal cells into cancer cells the immune systems reveal different response 
mechanisms which are described in more detail in [211 [22] . Firstly, the nascent transformed 
cells have to be identified. Candidates for the detection of tumor cells are the components of 
the innate immune system known as Natural Killer cells (NK), Natural Killer T-cells (NKT) 
and so called 7(5 T-cells. In case the tumor cells have been recognized the killer cells produce 
the cytokine Interferon-7 (IFN-7) as an important immunologic regulator [2^. Moreover, 
IFN-7 can cause the death of the tumor cell directly via apoptotic mechanisms [26]. The 
released IFN-7 leads to a stimulation of both the innate (activation of macrophages and 
presentation of antigens by dendrite cells (DC)) and the adaptive immune response (genera- 
tion of antigen-specific B- and T-lymphocytes). Eventually, the lymphocytes (CD8-positive 
T cells) migrate to the tumor site, detect the tumor cells and initiate a powerful immune 
reaction which may end up in the destruction of the tumor tissue. The complete suppres- 
sion of the cancer by the immune system is only one scenario. Likewise an imprinting of 
the tumor cells by their immunologic environment can occur during the tumor-immune cells 
reaction. So a selective pressure is exerted on the tumor which favors the creation of tumor 
cell clones that offer a low or even a non-immunogenic behavior. The very different response 
reflects the paradox role of the immune system of cancer promotion due to a sculpting of the 
immunogenic phenotype of the tumor. The numerous genetic alterations of the cancer cells 
during the sculpting process can be regarded as a sequence of stochastic events. Therefore, 
the modeling of the situation in a mathematical model should include both deterministic 
and stochastic parts. In addition the hypothesis of immunoediting suggest the occurrence 



of a phase with metastable states. Within this phase the tumor will neither grow to its final 
size nor be eliminated by the immune system. Because the tumor is under immunogenic 
control such a state can be regarded as tumor dormancy. As argued in [23] the period of 
this dormant state could be even of the life-time of the host. Despite of the short extract of 
possible effects one realizes that the immune system is a complex network where a variety 
of distinct cell types are involved with coordinated functions. An essential ingredient is that 
nearly all different cell types carry more than one functions. So the Natural Killer cells are 
able to release Interferon-7 and have simultaneously the ability to recognize and eliminate 
cancer cells. A further example is the immunomodulating agent IFN-7 which can on the 
one hand promote the proliferation of lymphocytes and on the other hand can directly effect 
the life of a cancer cell. 

Due to this diversity of cells and their functions and the fact that the interplay between 
tumor and immune cells is far from being understood completely the development of a 
mathematical model is necessary. Although one cannot expect that such models cover all the 
underlying biological aspects. Especially a very detailed description of the tumor-immune 
interaction seems not to be realistic. Otherwise such a coarsened model should include the 
main features of the immune system, namely detection, stimulation and elimination of tumor 
cells. Our approach simulates the different functions by introducing two kinds of immune 
cells named tumor detecting cells (TDC) and effector cells. The detecting cells are able 
to recognize the malignant cancer cells and additionally they stimulate the production of 
effector cells. The last ones have the ability to kill tumor cells. Insofar we map the three 
functions of the immune system onto two artificial cell types, the detecting cells and the 
effector cells. This mapping of the main functions of the immune cells allows us to construct 
a mathematical model the details of which are discussed in the following section. 

III. MODEL 

As discussed in the previous section, the immune system of the human body comprises 
various components which interact mutually. Moreover, the tumor cells are subjected to 
genetic alterations. Therefore, the tumor system can be regarded as being composed of 
different kinds of cells. In order to present an accessible theoretical model of a possible 
immune reaction against tumor growth we refer to the following coarsened description. The 



tumor system is assumed to consist of one single cell type the density of which at time t is 
denoted as X{t). Unlike the immune system is realized by two kinds of cells responsible for 
detection, stimulation and elimination, respectively. The elimination process is performed by 
the effector cells with density Y{t) which are able to kill the tumor cells. The second immune 
cell type -the tumor detecting cells (TDC) designated as Z{t)- have the ability to recognize 
the harmful cancer cells and in addition stimulate the proliferation of the effector cells. As 
basic model the three-state system obeys the following set of deterministic equations 

— X(t) =a (X(t) - bX^{t)) - cX{t) Y{t) , 
jY{t)=eY{t)Z{t)-~pY{t), (1) 

where the parameters a, 6, c, e, p, fl > will be discussed now. This model incorporate a 
logistic growth of the cancer cells X{t) with the birth rate a. The undisturbed evolution 
of the cancer would end when the tumor reaches its final size -the carrying capacity b~^. 
The effector cells Y{t) can interact with the tumor cells and hence the size of the tumor is 
reduced. The parameter c is a measure for the strength of the tumor-effector cell reaction. 
As suggested in the previous section the effector cells with the ability to kill the cancer 
cells do not exist without an external stimulus. The production of the effector cells will 
be mediated by the TDC with density Z{t). The term eY{t) Z{t) in Eq. ([I]) describes the 
initiation of effector cells due to the TDC. The parameter e is the production rate. Because 
the immune system can exert their influence only for a limited period, we have introduced 
the terms —pY{t) and —jlZ{t) in Eqs. ([I|. They reflect the finite lifetime p~^ and jl~^ of the 
effector cells and the TDC, respectively. As visible from Eq. ([l]) an elimination of the tumor 
is not possible within this approach because a release-term for the tumor-recognizing cells 
is not taken into account and thus effector cells are not produced. The three-state model 
for (X(t), F(t), Z{t)) in Eq. ([I]) offers two stationary states (0, 0, 0) and {h^^, 0, 0) where the 
tumor-free state X^ = is never stable. Instead of that the state X^ = X{t — )■ oo) = b^^ 
with a finite tumor population is realized. Eq. ([T| predicts that the tumor will always reach 
its final size determined by the carrying capacity. As discussed in Sec. II the tumor-immune 
interaction is subjected to numerous stochastic events. In the following we will demonstrate 
that random forces are able to create a birth term for the TDC Z{t). As the consequence 
the behavior of the system is changed drastically. To reduce the number of parameters let 



us introduce dimensionless variables according to 



x = hX^y='-Y^z='-Z^p=P-,,- , 
a a a a 



t = at . 



(2) 



In terms of these quantities and under introducing random forces i]i{t) the deterministic set 
of Eq. ([I]) is changed to the stochastic differential equations 



dt 



x,(t)=V'i[x(t)]+ai[x(t)]r/,-(t). 



(3) 



Here for simplicity of notation the dimensionless time variable i is replaced by t and sum- 
mation over double indices is understood. Eq. (Is]) describes the noisy tumor-immune inter- 
action. The vector i/) and the matrix f2 are defined by 
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Further, we have introduced the vector x = (x, y, 2;) and the vector of the stochastic force 
T] = {r]x,r]y,r]z), i.e. the noise r]i is associated with the cell type Xj. Eq. ^ and Eq. Q 
include the obvious possibility that the tumor cells x are coupled to the random force r]z 
originated in the TDC subsystem. This coupling term appears in the equation of motion of 
the TDC z, too. Because the tumor itself is thought to be a source of stochastic influences. 
So the couplings supposed between 2;-cells and the noise force 7]^ stemming from the tumor 
cells. Such a coupling term occurs in the evolution equation of the cancer cells x as well as 
in that one of the TDC z. The special form of the couplings was chosen to emphasize the 
importance of recognizing the tumor cells and the according stochastic events. The noisy 
properties are expressed by the following relations 
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The components rik{t) have a zero mean. In the limit that the correlation time tends to 
zero, Tkl -^ 0, the usual white noise properties are recovered. The correlation strength and 



correlation time matrices D and r, respectively, are assumed to take the forms 
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The matrix of the correlation time r^i reveals that all auto-correlations are characterized by 
the finite correlation time r whereas the cross-correlation functions with strengths R, S and 
P offers white noise properties with the (5-function according to Eq. (pi). 

IV. PROBABILITY DISTRIBUTION 



In this section we derive the probability distribution P(x, t) which is related to the set 
of stochastic equations determined by Eqs. (|3])-([6]). Following [271 128] we define 



P(x,t) = (,5[x(t)-x]>. 



(7) 



Here the < ... > means the average over all realizations of the stochastic process. The 
vector x(t) represents the stochastic process whereas the x are the possible realizations of 
the process at time t. Due to the colored noise the corresponding Fokker-Planck equation 
can be obtained only approximatively in lowest order of the correlation time. The time 
evolution of Eq. (JTl) can be written in the form 

d 



dt 



Pi^,t) =CP{x,t) 



(8) 



In deriving this expression we have used the time evolution of x(t) according to Eq. ([s]), the 
Novikov theorem |29] and the correlation function in Eq. ([s]) with Tki and Dki presented in 
Eq. ([6]). The form of the operator C is given in a correlation time and cumulant expansion 
[30H32] by 
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The single probability distribution is determined by the operator £ in Eq. (|9j). Notice that 
the representation is valid for sufficiently large times scale compared with the correlation 



times T^i when transient terms are negligible. Eqs. (pj)-(lO) enable to find the equation of 
motion for the expectation values {xj{t)). It follows 
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Remark that in the limiting case of white noise all terms including tm vanish. Further 



we want to point out that the expression in Eq. (11) contains quadratic terms like {xiXj) 
due to the nonlinear system in Eq. (|3|. In the same manner as before one can derive a 
higher order joint probability distribution, see [23] • Following this procedure we get a whole 
hierarchy of evolution equations. Instead of that let us make the simplest approximation 



(xiXj) = {xi){xj). Under this approximation Eq. (11) and by applying Eq. (10) the equation 
of motion for the mean values can be rewritten as 

|(x(t)) =[1 + Ril - D^t) + \d.,D,t] (x(t)) - {x{t)f - (x(t)) {y{t)) + [D,.(l + Rr)] {z{t)) 

j^ivit)) ={y{t)) {z{t)) - {p - Dy) {y{t)) , 

j^{z{t)) =[R{1 - D^t) + D^D^r] (x(t)) - [/i - (i?(l + D^r) + D.,{1 + ^D,r))] {z{t)) . 

(12) 



As can be seen from Eq. ( 12 ) the random process referring to the correlation strength and 
correlation time presented in Eq. dol) infiuences the dynamical system in a significant manner. 
The behavior is illustrated in Fig. [T] Let us compare the results of the stochastic approach 
with the deterministic model. The birth rate of the tumor cells (x) are affected by the noise 
correlation strengths D^, D^ associated with the tumor cells and the tumor detecting cells 
as well as their cross-correlation R. Likewise the decay rate p in the equation for the effector 
cells {y) is reduced by the noise strength Dy, i.e. by the noise related to the effector cell 
subsystem itself. As a new nontrivial result we find noise induced terms in the evolution 



equation Eq. (12). So there appears a term oc {x{t)) in the equation for the z-cells which are 



able to recognize cancer cells. In the same manner the generating term oc {z{t)) arises in the 
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FIG. 1. (Color online) Schematic illustration of the model presented in Eq. (12) 



equation for the tumor cells. These terms are originated exclusively due to the randomness 
in Eqs. (IS]) and Q. In this context the most important parameter is played by the cross 
correlation strength R of the correlation function Xxz = Xzx, i-e. the correlation between 
the noise sources inherent in the tumor cells and the tumor detecting cells. Notice that 
in the noise-free case such an interlink between these two cells is missing, compare again 
Fig. [T| Moreover, the death rate /i is altered due to the stochastic process. Based on the 



implementation of noisy forces the resulting Eq. ( 12 ) differs from the deterministic equation 



twice, (i) Firstly, the birth and death rates as well are altered due to stochastic parameters 
such as the correlation time r and the correlation strengths D^, Dy, Dz and R defined in 
Eq. ^. Although the cross-correlations S and P were included in Eqs. ^ and (|4]), they 



do not appear in the final expression Eq. (12). This fact is related to the special kind of 



multiplicative noise of our model, (ii) Secondly, two new terms exist in Eq. (12). The origin 



of both can be solely ascribed to stochastic sources. Regarding the evolution of the tumor 
detecting cells {z) the new term disappears in case R = and r = 0. So both parameters 
R -the cross-correlation strength between 7]^ and riz and r -the correlation time of the auto- 
correlation functions {i]i(t)Vi(t'))' ^^^ of significant relevance. In the subsequent section the 
analysis is focused in particular on both parameters. 
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V. RESULTS AND DISCUSSION 



As remarked in the previous section the parameters R and r inhere a special meaning in 



discussing the set of Eq. (12). Before we proceed with the stabihty analysis and the results 
we want to estimate the model parameters. The starting point is the deterministic Eq. ([I]). 
One finds different values for the intrinsic tumor growth rate a: 0.18 day~^ [Ij and 0.51 day~^ 
[1]. Our own study leads to 0.57 day~^ [12j. The first two values are based on mouse models 
while the latter one was obtained by means of in vitro cultivation of tumor cells. The growth 
rate is insofar of importance as it determines the time scale of the dynamics, see Eq. ^ 
{t [in days] = t/a). Here we use a = 0.5 day^ . Thus, t = 1 is tantamount to t = 2 days. 
An estimation of the carrying capacity is b^^ = 10^ cells [H H]. Further, the reaction rates 
take approximately c = 10~^ cell~^ day~^ = e [H [21 SI E]- An estimation for the decay rates 
in Eq. ([I]) is given by p = 3 x 10~^ day~^ and fl = 10 day~^ [21 II]- In relating our results 
to real units one should take into account the scaling properties Eq. pi). All the results 
are collected at the end of this section in Tab. [H For the subsequent analysis it is more 
convenient to use dimensionless quantities. The both most relevant parameters of stochastic 
forces are the auto-correlation time r and the cross-correlation strength R. Both quantities 
R and r will be altered within the interval [0,5]. The remaining parameters are assumed 
to be fixed, i.e. D^ = 2.1, D^ = 1.2 and Dy = 0.01. The values for D^ and D^ are chosen 
arbitrarily, whereas the value for Dy is suggested to be smaller than p = p/a = 0.06 in order 



to guarantee a sufficient stability of the differential equation system Eq. (12), cf. the term 
oc {y). Moreover, since we consider cell populations the solutions of Eq. ([I]) should yield 
positive values for the cell numbers. So values of R and r are excluded when they induce 
negative values for the cell populations. 

Now we perform the stability analysis according to the tumor-immune cells reaction 



system satisfying Eq. (12). We note that the numerical bifurcation analysis is performed by 
means of the program |34j which contains the bifurcation tool [35J. This set of equations 
exhibits three different equilibria, i.e. the tumor-free Ei = (0,0,0), and two non-tumor-free 
states designated as E2 and E3. The last ones are given by lengthy expressions in terms 
of the model parameters. Only one of the three equilibria is stable simultaneously. It is 
also possible that the total system becomes unstable as discussed below. The solution of 



Eq. (12) depends strongly on the correlation time r and the cross-correlation strength R. 
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FIG. 2. Behavior of the solution representing the regions I-III mentioned in the text. The param- 
eters take p = 0.06, /x = 20, Dx = 2.1, Dy = 0.01 and Dz = 1.2. (a)-(c): Bifurcation diagrams, 
(d) Limit cycles in the (x) — (z)— pahse space. 



Concerning r we find three different regions (labeled as I-III) where the solution of Eq. (12) 



has different properties. The threshold values referring to our specific numerical values of 
the remaining model parameters are r^^ = 0.636 and Tc^ = 4.016. We proceed by considering 
these three regions determined by the correlation time r G [0,5]. As fixed initial values for 
the tumor and the tumor detecting cells, respectively, we choose {x(t = 0)) = 10~^ and 
{z(t = 0)) = . This reflects a situation where the tumor is small and tumor detecting 
cells are not present. In our case this equals an initial tumor cell number of 10'^ cells which 
is clinically not detectable (early stages of tumor evolution). Tumor detecting cells should 
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be only generated due to stochastic influences and did not exist a priori. For the initial 
state of the effector cell population we have to choose a non-zero value of {y(t = 0)) > 0. 
Otherwise the solution will differ from the predictions based on the stability analysis which 



is due to the structure of the differential equation system in Eq. (12). In some regions this 
value can be very small (almost zero, e.g. ~ 10~^°) whereas in other areas corresponding to 
the parameters R and r the stability of the solution depends on {y{t = 0)) signiflcantly . In 
what follows these distinct cases will also be discussed while we want to restrict the possible 
values for {y(t = 0)) to the left-open interval (0, 10]. 

Region / (0 < r < Tc^ = 0.636): Within this range a stable tumor-free state is missing 
for R G [0,5]. All three flxed points exist in this region. The solution tends either to the 
steady states E2 or E^ one of which is asymptotically stable which depends on the cross- 
correlation strength R. For instance in case of r = 0.3 the equilibrium value of the tumor 



cell population, designated as {x)q as a function of R, is depicted in Fig. 2(a) As is visible 
for all < i? < 5 the solution will always reach the flxed point E2. Further, the equilibrium 
value (x)o decreases with increasing R. In region I the initial value of the effector cells 



can take arbitrary positive values in (0, 10] without changing the solution of Eq. (12). An 



exemplary dynamical solution is illustrated in Fig. 3(a) 



Region II [tc-^^ = 0.636 < t < Tc^ = 4.016): In this area the behavior is changed and one 



observes diagrams like that one shown in Fig. 2(b) for r = 2.0. The three flxed points -^1,2,3 



survive for R G [0,5], but the stability is changed. If we start at i? = and increase the 
cross-correlation strength R the behavior of the solution traverses four different regions. For 
< i? < 1.242 the steady state E2 is stable. At i? = 1.242 a transcritical bifurcation occurs 
where E2 is not stable anymore. The flxed point E^ becomes stable but only within the 
interval 1.242 < R < 1.325. At i? = 1.325 another transcritical bifurcation occurs namely 
the transition to the tumor-free state Ei which becomes stable while E^ loses its stability. 
Biologically such a transition is of great relevance because it manifests that the immune 
system is able to eliminate a growing tumor provided the tumor-immune cells reaction is 
assisted by a cross-correlation between stochastic events occurring in the tumor and in the 



tumor detecting cells subsystem. Notice that the sector A in Fig. 2(b) has to be excluded 
because the eigenvalues of Ei = (0,0,0) develop an imaginary part indicating the solution 
tends to Ei on a stable spiral. However, during the evolution towards the equilibrium value 
{{x{t)) — )■ {x)o, t — )► 00) the tumor cell population {x{t)) takes negative values. This happens 
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FIG. 3. Exemplary dynamic solutions according to the regions I-III mentioned in the text. The 
parameters take p = 0.06, /i = 20, Dx = 2.1, Dy = 0.01 and Dz = 1.2. 



for R > 1.733 which is indicated by the sector A in Fig. 2(b) In the area R < 1.733 there 



are no restrictions on the initial value for the effector cells < {y{0)) < 10. In Fig. 3(b) 



the time evolution of the tumor cell number {x{t)) is shown for different values of the cross- 
correlation strength R. Summarizing the result we observe in region II the occurrence of 
tumor escape as well as the possibility of tumor elimination depending on the value of the 
cross-correlation strength R. 

Region III {tc^ = 4.016 < r < 5): In this parameter range one observes a new behavior 
determined by the cross-correlation strength R and the initial value of the effector cells 



{y(0)). For the following discussion we refer to Fig. 2(c) Starting from R = and increasing 
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description in the text. 



this parameter the solution of Eq. (12) tends to the stable fixed point E2. A related solution 



is represented in Fig. |3(c) where it needs a rather long time until E2 is reached, about 
t ~ 19000. The fixed point E2 is realized on a stable spiral within the interval < i? < 0.888. 
The smaller the cross-correlation strength R is the shorter is the time scale to reach E2, for 
instance we need t = 1000 in case R = 0.01. When the critical value R = 0.888 is exceeded a 



periodic limit cycle evolves related to the occurrence of a Hopf bifurcation. In Fig. 2(c) the 



minimal and the maximal numbers of tumor cells are plotted within such a limit cycle. The 
numerical values range below and above the former stable equilibrium E2 which becomes now 
unstable. After the Hopf bifurcation the steady states Ei and E^ are no longer detectable. 



Further, Fig. 2(c) reveals that the the parameter range is limited in which such stable 



periodic oscillations emerge. The dashed line represents the boundary to sector B, where 
the total system bifurcates into an unstable state and the dynamical system is uncontrollable 
anymore. Thus the sector B will be excluded as a domain of accessible solutions within our 
tumor-immune model. Nevertheless, periodic orbits can be observed for 0.888 < R < 1.101 
and fixed correlation time r = 4.5. For varying values of R the periodic solutions are 



depicted in the two-dimensional (x) — (2)— phase space, see Fig. 2(d), The numbers shown 
above each orbit is the frequency of the oscillations between two maxima. With growing 
R from 0.9 to 1.0 the minimal and the maximal cell numbers increase for both (x) and (z) 
while the frequency decreases, i.e. the period of the oscillation is enlarged. This result is 
also valid for the effector cells (y) which are not shown here. Coming back to the infiuence of 
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the initial values of the effector cells {y{0)). As already mentioned they play a decisive role 
in the regime of periodic limit cycle solutions. For very low initial values the system loses 
its stability and the solution is not accessible biologically. The curve which separates stable 
periodic cycles from unstable solutions is displayed in Fig. |4j The question arises what does 
it means for the real tumor-immune cells interactions if stable periodic oscillations occur? 
In that case we argue that an intensive interaction between the cancer and the immune cells 
takes place at the beginning of the tumor growth as well as after a long period later. If 
the nascent transformed cells start to grow up the immune system is able to detect this 
harmful process and it responds. Such an immune attack reduces the tumor size without to 



delete it completely. Notice that our numerical estimation represented in Fig. 2(c) is also 
compatible with an elimination of the tumor cells because the lower branch is partly not 
distinguishable from zero in the interval 0.888 < R < 1.101. But the tumor starts anew 
to grow up signalizing the latent facility that the tumor evolution goes on. Otherwise, if 
the cancer growth is continued the immune system remains active and consequently it is 
still able to eliminate a large amount of tumor cells. So after a certain time one expect 
that a balance between tumor growth and the response of the immune system is evolved. 
From here one concludes that the tumor is under the control of the immune system and 
a so called tumor dormant state emerges. In the same manner the region with R < 0.888 
(before the Hopf bifurcations appears) may also be interlinked to the tumor dormant state. 
In that case the number of tumor cells is low compared to other parameter regimes, see 
Fig. [2] For instance the value R = 0.6 yields an equilibrium tumor cell population of 5% of 



the carrying capacity in Fig. 2(c) However, as a result of our computations, the size of the 
maximal number of tumor cells within one cycle can take large values, e.g. for R = 1.1 we 
find {x)max = 0.88. The tumor reaches 88% of its carrying capacity. The time-dependent 



solution for R = 1.0 is depicted in Fig. 3(d) within the interval < t < 500. Eventually a 



periodic cycle with {x)max = 0.54 will be reached after t ~ 12000. 

At the end of this section we want to convert some dimensionless quantities into quantities 
with real units. The results are particularized in Tab. |Tj 

Not commented yet is that the strengths in the correlation functions in Eq. d5| carry the 
unit day" after conversion to real units. This follows from the correlation function in real 
units, i.e. {rjirj'j) a^ oc {Dija)/{Tij/a) where the intrinsic tumor growth rate a = 0.5 day~^ 
and Dij and Tij are given in arbitrary units. Thus, the strengths occurring in the noise-noise 
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quantity 




arbit 


rary units 


real units 


time t 








2 days 


auto-correlation time r 








2 days 


cross-correlation strength R 








0.5 day-i 


number of tumor cells {x) 








10^ cells 


number of effector cells {y) 








5 X 10^ cells 


number of tumor detecting cells 


{z) 






5 X 10^ cells 


frequency (period) of | 


:;ycles 






0.029 


« 0.01 day-^ 


according to Fig. 


2(d) 






(^ 


» 34.5) 


(ps 69 days) 



TABLE I. Comparison of model quantities in arbitrary and real units. See also Eq. Q. 
correlation functions in Eq. dsl) have the meaning of a rate. 

VI. CONCLUSIONS 

We have presented a mathematical model for the tumor-immune cells reactions which 
is essentially supplemented by stochastic forces. The parameters of the noise correlation 
function have a great impact on the behavior of the coupled tumor-immune cell interaction, 
especially on the response of the immune system. In particular we have emphasized that the 
auto-correlation time r and the cross-correlation strength R are able to control the evolution 
of the tumor. More precise, these two quantities discriminate whether the system tends to 
tumor suppression, tumor progression or tumor dormancy. The assistance of an inevitable 
noisy influences seems to play a crucial role during cancer genesis and growth in humans. 
The involved random forces may be originated within the tumor as well as inside the immune 
system and can even interact mutually which is manifested in the cross-correlation. Our 
model should be considered as an attempt toward a more detailed analysis of tumor-immune 
systems. But also the model studied elucidates that noise plays an decisive role in such 
systems. The model can be refined immediately, e.g. a finite correlation time is attributed 
to the cross-correlation functions, too. In that case the correlation time matrix, Eq. (pi), is 



modified and new terms in Eq.(12) occur. We believe that our approach includes the most 
relevant degrees of freedom. 
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